#include "KSADDetectorConstruction.hh"
#include "KSADGlobals.hh"
#include "G4Element.hh"
#include "G4Material.hh"
#include "G4Box.hh"
#include "G4Tubs.hh"
#include "G4LogicalVolume.hh"
#include "G4ThreeVector.hh"
#include "G4PVPlacement.hh"
#include <cmath>
#include <sstream>
#include "KSADRearDetector.hh"
#include "G4SDManager.hh"
#include "G4UImanager.hh"
#include "G4ios.hh"
#include "G4NistManager.hh"
#include "G4MaterialPropertiesTable.hh"
#include "G4OpticalSurface.hh"
#include "G4LogicalSkinSurface.hh"
#include "G4LogicalBorderSurface.hh"

KSADDetectorConstruction::KSADDetectorConstruction()
{
	;
}

KSADDetectorConstruction::~KSADDetectorConstruction(){;}

G4VPhysicalVolume *KSADDetectorConstruction::Construct()
{
	// Define materials
	
	G4NistManager* nistMan = G4NistManager::Instance();

	G4Element* N = new G4Element("Nitrogen","N",7,14.01*g/mole);
	G4Element* O = new G4Element("Oxygen","O",8,16.00*g/mole);
	G4Material* Air = new G4Material("Air",15*g/6.022/pow(10,23)/cm3,2);
	Air->AddElement(N,70*perCent);
	Air->AddElement(O,30*perCent);

	G4Material* CsI = nistMan->FindOrBuildMaterial("G4_CESIUM_IODIDE");
	G4MaterialPropertiesTable* CsIMPT = new G4MaterialPropertiesTable();
	
	const int scintEntries = 2;
	double energy[scintEntries] = {2*eV,4*eV};
	double rindex[scintEntries] = {1.82,2.08};
	double absorp[scintEntries] = {2*m,1*m};
	double scintFast[scintEntries] = {1,1};
	double scintSlow[scintEntries] = {2,2};

	CsIMPT->AddProperty("RINDEX",energy,rindex,scintEntries);
	CsIMPT->AddProperty("ABSLENGTH",energy,absorp,scintEntries);
	CsIMPT->AddProperty("FASTCOMPONENT",energy,scintFast,scintEntries);
	CsIMPT->AddProperty("SLOWCOMPONENT",energy,scintSlow,scintEntries);

	CsIMPT->AddConstProperty("SCINTILLATIONYIELD",100./MeV);
	CsIMPT->AddConstProperty("RESOLUTIONSCALE",1.0);
	CsIMPT->AddConstProperty("FASTTIMECONSTANT",1.*ns);
	CsIMPT->AddConstProperty("SLOWTIMECONSTANT",10.*ns);
	CsIMPT->AddConstProperty("YIELDRATIO",0.8);

	CsI->SetMaterialPropertiesTable(CsIMPT);

	G4Material* Al = nistMan->FindOrBuildMaterial("G4_Al");

	// Define world
	
	G4Box* worldSolid = new G4Box("worldSolid",2.5*m,2.5*m,5*m);
	G4LogicalVolume* worldLog = new G4LogicalVolume(worldSolid,Air,"worldLog");
	G4VPhysicalVolume* worldPhys = new G4PVPlacement(0,G4ThreeVector(0,0,0),worldLog,"world",0,0,0);

	// Define crystal block
	
	G4Tubs* crystalBlockSolid = new G4Tubs("crystalBlockSolid",0,0.95*m,25*cm,0,360*deg);
	G4LogicalVolume* crystalBlockLog = new G4LogicalVolume(crystalBlockSolid,Air,"crystalBlockLog");
	G4VPhysicalVolume* crystalBlockPhys = new G4PVPlacement(0,G4ThreeVector(0,0,25*cm),crystalBlockLog,"crystalBlock",worldLog,0,0);
	(void)crystalBlockPhys;

	// Define small crystals

	G4SDManager* SDMan = G4SDManager::GetSDMpointer();
	G4VSensitiveDetector* rearSensitiveDetector = new KSADRearDetector("rearSensitiveDetector");
	SDMan->AddNewDetector(rearSensitiveDetector);

	G4Box* smallCrystalSolid = new G4Box("smallCrystalSolid",crystalSize/2,crystalSize/2,25*cm);
	G4LogicalVolume* smallCrystalLog = new G4LogicalVolume(smallCrystalSolid,CsI,"smallCrystalLog");

	G4OpticalSurface* separatorOptSurface = new G4OpticalSurface("separatorOptSurface",glisur,polished,dielectric_metal);

	G4Box* smallCapSolid = new G4Box("smallCapSolid",crystalSize/2,crystalSize/2,0.025*mm);
	G4LogicalVolume* smallCapLog = new G4LogicalVolume(smallCapSolid,Al,"smallCapLog");
	G4LogicalSkinSurface* smallCapSurface = new G4LogicalSkinSurface("smallCapSurface",smallCapLog,separatorOptSurface);
	(void)smallCapSurface;
	G4VPhysicalVolume* smallCapPhys = new G4PVPlacement(0,G4ThreeVector(0,0,-25*cm+0.025*mm),smallCapLog,"smallCap",smallCrystalLog,0,0);
	(void)smallCapPhys;

	G4LogicalVolume* smallRearSensitiveLog = new G4LogicalVolume(smallCapSolid,CsI,"smallRearSensitiveLog",0,rearSensitiveDetector);
	G4VPhysicalVolume* smallRearSensitivePhys = new G4PVPlacement(0,G4ThreeVector(0,0,25*cm-0.025*mm),smallRearSensitiveLog,"smallRearSensitive",smallCrystalLog,0,0);
	(void)smallRearSensitivePhys;

	double RBPEnergy[2] = {2*eV,4*eV}; double RBPRefl[2] = {0.0,0.0}; double RBPEff[2] = {1.0,1.0};
	G4MaterialPropertiesTable* RBPMPT = new G4MaterialPropertiesTable();
	RBPMPT->AddProperty("REFLECTIVITY",RBPEnergy,RBPRefl,2);
	RBPMPT->AddProperty("EFFICIENCY",RBPEnergy,RBPEff,2);
	G4OpticalSurface* RBPOptSurface = new G4OpticalSurface("rearBackingPlateOptSurface",glisur,polished,dielectric_metal);
	RBPOptSurface->SetMaterialPropertiesTable(RBPMPT);

	G4Box* smallRBPSolid = new G4Box("smallRearBackingPlateSolid",crystalSize/2,crystalSize/2,0.0125*mm);
	G4LogicalVolume* smallRBPLog = new G4LogicalVolume(smallRBPSolid,Al,"smallRearBackingPlateLog");
	G4LogicalSkinSurface* smallRBPSurface = new G4LogicalSkinSurface("smallRearBackingPlateSurface",smallRBPLog,RBPOptSurface);
	(void)smallRBPSurface;
	G4VPhysicalVolume* smallRBPPhys = new G4PVPlacement(0,G4ThreeVector(0,0,0.0125*mm),smallRBPLog,"smallRearBackingPlate",smallRearSensitiveLog,0,0);
	(void)smallRBPPhys;
	
	G4Box* shortHorizSepSolid = new G4Box("shortHorizSepSolid",crystalSize/2,0.025*mm,25*cm-0.05*mm+0.0125*mm);
	G4Box* shortVertSepSolid = new G4Box("shortVertSepSolid",0.025*mm,crystalSize/2-0.05*mm,25*cm-0.05*mm+0.0125*mm);

	G4LogicalVolume* shortLeftSepLog = new G4LogicalVolume(shortVertSepSolid,Al,"shortLeftSepLog");
	G4LogicalSkinSurface* shortLeftSepSurface = new G4LogicalSkinSurface("shortLeftSepSurface",shortLeftSepLog,separatorOptSurface);
	(void)shortLeftSepSurface;
	G4VPhysicalVolume* shortLeftSepPhys = new G4PVPlacement(0,G4ThreeVector(-crystalSize/2+0.025*mm,0,0.0125*mm),shortLeftSepLog,"shortLeftSep",smallCrystalLog,0,0);
	(void)shortLeftSepPhys;

	G4LogicalVolume* shortRightSepLog = new G4LogicalVolume(shortVertSepSolid,Al,"shortRightSepLog");
	G4LogicalSkinSurface* shortRightSepSurface = new G4LogicalSkinSurface("shortRightSepSurface",shortRightSepLog,separatorOptSurface);
	(void)shortRightSepSurface;
	G4VPhysicalVolume* shortRightSepPhys = new G4PVPlacement(0,G4ThreeVector(crystalSize/2-0.025*mm,0,0.0125*mm),shortRightSepLog,"shortRightSep",smallCrystalLog,0,0);
	(void)shortRightSepPhys;

	G4LogicalVolume* shortTopSepLog = new G4LogicalVolume(shortHorizSepSolid,Al,"shortTopSepLog");
	G4LogicalSkinSurface* shortTopSepSurface = new G4LogicalSkinSurface("shortTopSepSurface",shortTopSepLog,separatorOptSurface);
	(void)shortTopSepSurface;
	G4VPhysicalVolume* shortTopSepPhys = new G4PVPlacement(0,G4ThreeVector(0,crystalSize/2-0.025*mm,0.0125*mm),shortTopSepLog,"shortTopSep",smallCrystalLog,0,0);
	(void)shortTopSepPhys;

	G4LogicalVolume* shortBottomSepLog = new G4LogicalVolume(shortHorizSepSolid,Al,"shortBottomSepLog");
	G4LogicalSkinSurface* shortBottomSepSurface = new G4LogicalSkinSurface("shortBottomSepSurface",shortBottomSepLog,separatorOptSurface);
	(void)shortBottomSepSurface;
	G4VPhysicalVolume* shortBottomSepPhys = new G4PVPlacement(0,G4ThreeVector(0,-crystalSize/2+0.025*mm,0.0125*mm),shortBottomSepLog,"shortBottomSep",smallCrystalLog,0,0);
	(void)shortBottomSepPhys;

	int smallCrystalDimension = 1.2*m/crystalSize;
	for(int i = 0;i < smallCrystalDimension;i++)
	{
		for(int j = 0;j < smallCrystalDimension;j++)
		{
			double xPos = (i-smallCrystalDimension/2+.5)*crystalSize;
			double yPos = (j-smallCrystalDimension/2+.5)*crystalSize;
			if(fabs(xPos) > 10*cm || fabs(yPos) > 10*cm)
			{
				std::stringstream ss;
				ss << "smallCrystal" << i << "," << j;
				G4String name = ss.str();
				G4VPhysicalVolume* smallCrystalPhys = new G4PVPlacement(0,G4ThreeVector(xPos,yPos,0),smallCrystalLog,name,crystalBlockLog,0,0);
				(void)smallCrystalPhys;
			}
		}
	}

	// Define large crystals

	G4Box* largeCrystalSolid = new G4Box("largeCrystalSolid",crystalSize,crystalSize,25*cm);
	G4LogicalVolume* largeCrystalLog = new G4LogicalVolume(largeCrystalSolid,CsI,"largeCrystalLog");

	G4Box* largeCapSolid = new G4Box("largeCapSolid",crystalSize,crystalSize,0.025*mm);
	G4LogicalVolume* largeCapLog = new G4LogicalVolume(largeCapSolid,Al,"largeCapLog");
	G4LogicalSkinSurface* largeCapSurface = new G4LogicalSkinSurface("largeCapSurface",largeCapLog,separatorOptSurface);
	(void)largeCapSurface;
	G4VPhysicalVolume* largeCapPhys = new G4PVPlacement(0,G4ThreeVector(0,0,-25*cm+0.025*mm),largeCapLog,"largeCap",largeCrystalLog,0,0);
	(void)largeCapPhys;

	G4LogicalVolume* largeRearSensitiveLog = new G4LogicalVolume(largeCapSolid,CsI,"largeRearSensitiveLog",0,rearSensitiveDetector);
	G4VPhysicalVolume* largeRearSensitivePhys = new G4PVPlacement(0,G4ThreeVector(0,0,25*cm-0.025*mm),largeRearSensitiveLog,"largeRearSensitive",largeCrystalLog,0,0);
	(void)largeRearSensitivePhys;

	G4Box* largeRBPSolid = new G4Box("largeRearBackingPlateSolid",crystalSize,crystalSize,0.0125*mm);
	G4LogicalVolume* largeRBPLog = new G4LogicalVolume(largeRBPSolid,Al,"largeRearBackingPlateLog");
	G4LogicalSkinSurface* largeRBPSurface = new G4LogicalSkinSurface("largeRearBackingPlateSurface",largeRBPLog,RBPOptSurface);
	(void)largeRBPSurface;
	G4VPhysicalVolume* largeRBPPhys = new G4PVPlacement(0,G4ThreeVector(0,0,0.0125*mm),largeRBPLog,"largeRearBackingPlate",largeRearSensitiveLog,0,0);
	(void)largeRBPPhys;
	
	G4Box* longHorizSepSolid = new G4Box("longHorizSepSolid",crystalSize,0.025*mm,25*cm-0.05*mm+0.0125*mm);
	G4Box* longVertSepSolid = new G4Box("longVertSepSolid",0.025*mm,crystalSize-0.05*mm,25*cm-0.05*mm+0.0125*mm);

	G4LogicalVolume* longLeftSepLog = new G4LogicalVolume(longVertSepSolid,Al,"longLeftSepLog");
	G4LogicalSkinSurface* longLeftSepSurface = new G4LogicalSkinSurface("longLeftSepSurface",longLeftSepLog,separatorOptSurface);
	(void)longLeftSepSurface;
	G4VPhysicalVolume* longLeftSepPhys = new G4PVPlacement(0,G4ThreeVector(-crystalSize+0.025*mm,0,0.0125*mm),longLeftSepLog,"longLeftSep",largeCrystalLog,0,0);
	(void)longLeftSepPhys;

	G4LogicalVolume* longRightSepLog = new G4LogicalVolume(longVertSepSolid,Al,"longRightSepLog");
	G4LogicalSkinSurface* longRightSepSurface = new G4LogicalSkinSurface("longRightSepSurface",longRightSepLog,separatorOptSurface);
	(void)longRightSepSurface;
	G4VPhysicalVolume* longRightSepPhys = new G4PVPlacement(0,G4ThreeVector(crystalSize-0.025*mm,0,0.0125*mm),longRightSepLog,"longRightSep",largeCrystalLog,0,0);
	(void)longRightSepPhys;

	G4LogicalVolume* longTopSepLog = new G4LogicalVolume(longHorizSepSolid,Al,"longTopSepLog");
	G4LogicalSkinSurface* longTopSepSurface = new G4LogicalSkinSurface("longTopSepSurface",longTopSepLog,separatorOptSurface);
	(void)longTopSepSurface;
	G4VPhysicalVolume* longTopSepPhys = new G4PVPlacement(0,G4ThreeVector(0,crystalSize-0.025*mm,0.0125*mm),longTopSepLog,"longTopSep",largeCrystalLog,0,0);
	(void)longTopSepPhys;

	G4LogicalVolume* longBottomSepLog = new G4LogicalVolume(longHorizSepSolid,Al,"longBottomSepLog");
	G4LogicalSkinSurface* longBottomSepSurface = new G4LogicalSkinSurface("longBottomSepSurface",longBottomSepLog,separatorOptSurface);
	(void)longBottomSepSurface;
	G4VPhysicalVolume* longBottomSepPhys = new G4PVPlacement(0,G4ThreeVector(0,-crystalSize+0.025*mm,0.0125*mm),longBottomSepLog,"longBottomSep",largeCrystalLog,0,0);
	(void)longBottomSepPhys;

	int largeCrystalDimension = 1.9*m/crystalSize/2;
	for(int i = 0;i < largeCrystalDimension;i++)
	{
		for(int j = 0;j < largeCrystalDimension;j++)
		{
			double xPos = (i-largeCrystalDimension/2+.5)*2*crystalSize;
			double yPos = (j-largeCrystalDimension/2+.5)*2*crystalSize;
			if(fabs(xPos) > 0.6*m || fabs(yPos) > 0.6*m)
			{
				if(sqrt(pow(fabs(xPos)+crystalSize,2)+pow(fabs(yPos)+crystalSize,2)) < 0.95*m)
				{
					std::stringstream ss;
					ss << "largeCrystal" << i << "," << j;
					G4String name = ss.str();
					G4VPhysicalVolume* largeCrystalPhys = new G4PVPlacement(0,G4ThreeVector(xPos,yPos,0),largeCrystalLog,name,crystalBlockLog,0,0);
					(void)largeCrystalPhys;
				}
			}
		}
	}

	return worldPhys;
}
